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ABSTRACT 

The unstable invariant manifolds of the short-period family of periodic orbits around the unstable Lagrangian points Lj and L2 of a 
barred galaxy define loci in the configuration space which take the form of a trailing spiral pattern. In previous works we have explored 
the association of such a pattern to the observed spiral pattern in A'-body models of barred-spiral galaxies and found it to be quite 
relevant. Our aims in the present paper are: a) to investigate this association in the case of the self-consistent models of Kaufmann & 
Contopoulos (1996) which provide an approximation of real barred-spiral galaxies, b) to examine the dynamical role played by each 
of the non-axisymmetric components of the potential, i.e. the bar and the spiral perturbation, and their consequences on the form of 
the invariant manifolds, and c) to examine the relation of 'response' models of barred-spiral galaxies with the theory of the invariant 
manifolds. Our method relies on calculating the invariant manifolds for values of the Jacobi constant close to its value for Li and 
L2. Our main results are the following: a) The invariant manifolds yield the correct form of the imposed spiral pattern provided that 
their calculation is done with the spiral potential term turned on. We provide a theoretical model explaining the form of the invariant 
manifolds that supports the spiral structure. The azimuthal displacement of the Lagrangian points with respect to the bar's major axis 
is a crucial parameter in this modeling. When this is taken into account, the manifolds necessarily develop in a spiral-like domain 
of the configuration space, delimited from below by the boundary of a banana-like non-permitted domain, and from above either by 
rotational KAM tori or by cantori forming a stickiness zone. On the contrary, if the whole non-axisymmetric perturbation is artificially 
'aligned' with the bar (i.e. there is no azimuthal shift of the Lagrangian manifolds), the manifolds support a ring rather than a spiral 
structure, b) We construct 'spiral response' models on the basis of the theory of the invariant manifolds and examine the connection 
of the latter to the 'response' models (Patsis 2006) used to fit real barred-spiral galaxies, explaining how are the manifolds related to 
a number of morphological features seen in such models. 
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1. Introduction 

The ordered or chaotic nature of orbits in barred galaxies 
has been the subject of many investigations in the literature 
(Contopoulos 1981; Pfenniger 1984; Sparke & Sellwood 1987; 
Pfenniger & Frendli 1991; Kaufmann & Contopoulos 1996; 
Patsis et al. 1997; Fux 2001; Pichardo et al. 2004; Kaufmann 
& Patsis 2005). Interest to this problem stems from the fact that 
the existence (and degree) of chaos has direct consequences on 
the morphological features of a rotating galaxy. In particular, the 
appearance of a large degree of chaos in the corotation region is 
one of the main reasons for why the bars terminate near corota- 
tion (Contopoulos 198 1, see Contopoulos 2002, pp. 473-474 for 
a review). 

Beyond corotation, prominent structures such as rings or spi- 
ral arms are commonly observed. The role of the chaotic orbits 
in the dynamics of such structures is still a widely open prob- 
lem, but recently some progress was made towards its under- 
standing. In particular, a theoretical model has been proposed 
and numerically explored (Voglis et al. 2006a, 2006b; Romero- 
Gomez et al. 2006, 2007), according to which the spiral arms (or 
rings) are supported by the unstable invariant manifolds of the 
two short period families of unstable periodic orbits around the 
unstable Lagrangian equilibria Li and L2 (called hereafter the 
PLl and PL2 families respectively). This theory was extended 



by Tsoutsis et al. (2008), by examining the contribution of the 
unstable manifolds of other families, besides PLl or PL2, to the 
same phenomenon. The importance of the chaotic orbits in sup- 
porting the spiral structure of barred galaxies has also been em- 
phasized by Patsis (2006). 

The following is a brief account of the theory of the invariant 
manifolds: 

1) We consider a 2D approximation of the orbits in the disk 
plane of a barred - spiral galaxy, given by the Hamiltonian 

1 

H{r,e,pr,pe)=-{p^r + -;j)-^pPe + Vo{r) + Vi{r,e) = Ej .(1) 

In this expression, (r,ff) are polar coordinates in the rotating 
frame, pr - r, pg - r^(6 + Q.p) is the angular momentum in 
the rest frame, Vo is the axisymmetric potential and Vi is the 
non-axisymmetric potential perturbation due to the bar and to 
the spiral arms. O,, is the angular speed of the rotating frame, 
which coincides with the bar-spiral pattern speed in an approxi- 
mation in which the latter is assumed to be unique. 

2) The Hamiltonian flow under ([TJ yields two stable (L4, Lg) 
and two unstable (Li, L2) Lagrangian equilibrium points (in the 
rotating frame) at which a star corotates with the pattern. The 
unstable manifold "W^^ of L[ is defined as the set of all the ini- 
tial conditions (ro, Oo,pro,pBo) in the phase space for which the 
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resulting orbit tends asymptotically to Li in the backward sense 
of time, namely 

lim ||e(f; ro, 0o, /?rO, - (^l, , 0l, , 0, = O) (2) 

where Qit, rQ,OQ, pro, peo) denotes the position (point in phase 
space) at time f of a particle along an orbit starting with 
the above initial conditions, and the norm || ■ || means the 
Euclidean distance between this point and the phase space point 
(ri, , , 0, Q.p), corresponding to Li . All the points of the man- 
ifold yield the same value of the Jacobi constant, equal to 
Ej - Ej^Li. Furthermore, since Li is simply unstable, "W^^ is 
a two-dimensional manifold embedded in the three-dimensional 
hypersurface of the phase space corresponding to a fixed Jacobi 
constant Ej = E^u. Similar definitions and properties hold for 
L2 and "W^^ and for the stable manifolds , i.e. the sets 
of initial conditions tending asymptotically to Li, or L2 in the 
forward sense of time, as f — > 00. 

3) For Ej > EjLi, a short-period unstable periodic orbit 
(PLl) bifurcates from Li (and the symmetric orbit PL2 from 
L2). This orbit forms a small loop around Li (Fig. la, thick solid 
curve), which corresponds to a ID-torus in the phase space. This 
torus is 'whiskered', i.e., it possesses its own asymptotic mani- 
folds. In particular, the unstable manifold of PLl is now defined 
as 

^PLi = [[j(fo, 00, Pro, Peo) : 

lim \\Q(t;ro,6o,p,o,pgo)-PLl\\ = o! (3) 

where the notation || - || refers to the minimum of the distances 
of Q(t) from the locus of all the phase space points of the or- 
bit PLl. For any fixed value of £7 > Eju, "W^^^ is a two- 
dimensional manifold embedded in the three-dimensional hy- 
persurface of constant Ej. Figure la shows the projection of a 
small part of this manifold, close to PLl, in the configuration 
space X - r cos 0,y - r sin Q. This is drawn approximately, by 
calculating a number of orbits with initial conditions on "Wp^^ , 
and close to PLl. The possibility to find such initial conditions 
is guaranteed by the fact that the manifold ^^^j is tangent to 
the unstable manifold of the linearized Hamiltonian flow near 
PLl (the so-called Grobman 1959 and Hartman 1960 theorem), 
and the latter is calculated by diagonalizing the Floquet matrix 
of the orbit PLl. In Fig. la we draw the part of the manifold ly- 
ing outside corotation for a particular model of barred galaxy. 
We can see that close to PLl the orbits form epicyclic loops of 
size nearly equal to the PLl loop, while, in the same time, the 
guiding center recedes from PLl along a path which yields a 
trailing spiral arm. An analysis of the linearized flow yields that 
the deviation of the guiding center from PLl is exponential in 
time, with a rate determined by the positive characteristic ex- 
ponents of the Floquet matrix of PLl. Furthermore, in generic 
galactic potentials all the orbits on 'Wp^^ are chaotic. (The same 
phenomena hold for the orbit PL2 and the manifold "Wp^^ also 
plotted in Fig. la) 

4) In strongly nonlinear models (as is the case of strongly 
barred galaxies with conspicuous spiral arms), further away from 
PLl the size of the epicycles becomes great (it may exceed the 
size of the bar). Such an example is shown in Fig. lb, referring 
to the orbits of the "Wp^j family in a A^-Body model of a barred 
galaxy (Voglis et al. 2006a, 2006b). We see that one such orbit 
(bold) forms two relatively small loops near PLl, reaching the 



apocentric positions Ai and A2, but the exponential recession of 
the guiding center is so fast that there is no loop formed between 
the second and third (A3) apocentric positions. Furthermore, the 
fourth apocentric position is at a distance about twice the bar's 
major semi-axis. Further integration beyond that of Fig. lb shows 
that, in fact, all these orbits belong to the so-called 'hot popula- 
tion' (Sparke & Sellwood 1987), i.e., the orbits make several 
consecutive oscillations in and out of corotation. Kaufmann and 
Contopoulos (1996, their figure 21a) suggested that such orbits 
can partly support the bar and partly the spiral arms. 

To understand how the chaotic orbits may establish a long- 
time flow supporting the spiral structure, Voglis et al. (2006a) 
examined a particular subset of points of the unstable mani- 
fold "Wp^^ , namely the locus of all apocentric positions of the 
orbits on '^p^- This is defined by taking the intersection of 
'Wp^j with the so-called surface of section of the apocentric po- 
sitions on a hypersurface of constant Ej, defined as the set of 
phase space points satisfying pr = 0, p,- < 0. We stress that the 
choice of such a surface of section proves to be very relevant to 
the particular type of study undertaken here. Indeed, whenever 
Pr - 0, Eq.([T]) yields points {r,9,pg) which can be projected 
either on the {9, pg) plane (called the 'phase portrait' by virtue 
of the fact that the variables (6, pg) are canonically conjugated), 
or the plane (0, r), i.e. the usual plane of motion in the rotating 
frame. However, in typical galactic potentials the variables r and 
pg at the apocenters have a monotonic relation, and this implies 
that the phase portraits and (r, 0) portraits for this particular type 
of surface of section are isomorphic. This allows us to unravel 
immediately the consequences of the phase space dynamical fea- 
tures, as seen in the phase portraits, to the morphological features 
of the system, as seen in the usual disk plane of motion with co- 
ordinates {x,y) - (rcos 0,rsm0). Such a comparison would not 
be possible by the use of a traditional surface of section, such as 
(x, x) for y = 0, j > (or j < 0). 

Returning to the role of the invariant manifolds, the inter- 
section of ^p£j with the apocentric surface of section yields 
an one-dimensional locus of points. Such a locus can be pro- 
jected on either the phase portrait plane (0, pg) or the configura- 
tion space (r, 0). Figure Ic shows the latter projection in the case 
of the same manifold as in Fig. lb, but calculated for a much 
larger length. Every point in Fig.lc corresponds to one apoc- 
entric position of a chaotic orbit with initial conditions on the 
unstable manifold. Clearly, the apocentric positions along "Wp^j 
yield a locus which also supports a trailing spiral arm over, how- 
ever, a much larger extent of 'Wp^^ than in the case of Fig. lb. 
The manifold of Fig. Ic takes a typical form known in dynamical 
systems' theory to be associated with the so-called phenomenon 
of homoclinic chaos. Briefly, the manifold develops lobes form- 
ing oscillations close to the apocentric points of the periodic or- 
bits PLl or PL2. Such oscillations are analyzed in detail in the 
sequel. 

We should stress that an analysis of the Floquet matrix of 
the PLl or PL2 families yields that only the directions of the 
unstable invariant manifolds 'Wp^j , '^pi2 s^^ch as to define 
trailing spiral arms, while, close to Li or L2, the stable mani- 
folds 'Wp^^, "V^^pui define leading spiral arms. Furthermore, in 
the forward sense of time the chaotic orbits are attracted in di- 
rections of the phase space along the unstable manifolds. In the 
sequel we no longer refer to the stable manifolds 'Wp^^, 'Wpi2, 
and the term 'invariant manifolds' always implies the unstable 
manifolds . 

In summary, the theory of the invariant manifolds, viewed as 
either the loci on which lies the continuous flow of a swarm of 
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Fig. 1. (a) Projections of the invariant manifolds "Wp^^ and 'Wp^2 the configuration space, approximated by an ensemble of 20 
orbits with initial conditions on the unstable direction of the tangent space to the periodic orbits PLl and PL2 (small bold circles) 
in the neighborhood of Li and Lj. The potential used is from model A (see Sect. 2). The bar is aligned to the y-axis and rotates 
counterclockwise, (b) Same as in (a) but in a strongly non-linear model (in this case we use the potential of the A^-Body simulation 
analyzed in Voglis et al. 2006a. Only the manifold ^p^i plotted), (c) The projection of the intersection of the manifolds "Wp^^ 
and "Wp^j ^^^^ apocentric surface of section pr = 0, pr <0 on the configuration space of the model of (b). 



orbits (Romero-Gomez et al. 2006, 2007), or the loci of apoc- 
entric positions of these orbits (Voghs et al. 2006a; Tsoutsis et 
al. 2008), predicts the formation by the manifolds of a trailing 
spiral pattern beyond corotation. Naturally, the central question 
that should be posed now is whether (and up to what extent) the 
spiral arms formed self-consistently in real galaxies can be asso- 
ciated with the spiral patterns formed by the invariant manifolds 
"^PLi 2- previous works (Voglis et al. 2006a; Tsoutsis et 

al. 2008), we examined this question by considering the spiral 
arms formed in an A^-Body model of a barred galaxy and found 
such an association to be quite relevant. 

In the present paper, our main goal is to examine the same 
question in simple models of real barred-spiral galaxies for 
which some reliable estimation of both the gravitational poten- 
tial and the pattern speed have been provided in the literature 
by methods independent of the previous considerations. To this 
end, we selected the potential models and pattern speeds re- 
ported in the study of Kaufmann & Contopoulos ( 1 996) for three 
real galaxies, NGC3992, NGC1073 and NGC1398. This choice 
is motivated by the fact that Kaufmann & Contopoulos (1996) 
constructed approximate self-consistent models of the studied 
galaxies based on the response density of the superposition of 
many stellar dynamical orbits. Thus, their study yielded not only 
plausible values of the potential parameters, or the pattern speed, 
but also the decomposition of the potential into components, i.e., 
Vhaio, Vdisk, Vhar and Vspirai- This allows us to check the role of 
each of these components, in particular of the non-axisymmetric 
ones Vbar and Vspirai, in the theory. It should be noted that the 
self-consistent technique, pioneered by Schwarzschild (1979), 
has been used extensively to provide reliable models of galax- 
ies, despite the fact that there is no a priori guarantee of the sta- 
biUty of such models that should ideally be probed via A^-body 
simulations (see e.g. Smith and Miller 1982). 

Besides re-confirming that the invariant manifolds do corre- 
late well with the spiral arms found in the self-consistent mod- 
els of Kaufmann & Contopoulos (1996), our investigation led 
to a second non-trivial result analyzed in detail in the sequel: 
In all three models the bar component is dominant over the spi- 
ral component within a large radial extent, but not in a narrow 
zone beyond corotation. This implies that if one uses only the 



bar component to calculate the manifolds, the latter yield ring 
rather than spiral structures. Furthermore, if one adds the spi- 
ral perturbation to the potential, but gives no azimuthal tilting 
to the associated m = 2 Fourier component, the manifolds be- 
come more open as regards their radial extent, but remain quite 
symmetric as regards their orientation with respect to the bar's 
major axis, thus still defining rings rather than spiral arms. Only 
when the azimuthal deformation of the equipotential surfaces 
due to a really spiral-Uke perturbation is taken into account (in 
the Kaufmann & Contopoulos (1996) paper this was modeled as 
a simple logarithmic spiral), the manifolds are found to follow 
closely the spiral arms of the self-consistent models. In some 
numerical experiments (see Sect. 3 below) we managed to ob- 
tain a kind of spiral pattern formed by the initial segments of 
the invariant manifolds in pure bar models, having, however, to 
drastically depart from the bar parameters given in Kaufmann 
and Contopoulos' self-consistent models, and pushing the bar's 
amplitude to highly non-physical values. But even in that case, 
the manifold-induced spiral arms are quite difi'erent from the spi- 
ral arms of the self-consistent models, and they disappear when 
the manifolds are computed for a longer length. Such an inves- 
tigation demonstrates that while in principle the strength of the 
quadrupole moment of the bar's potential causes a 'thickening' 
of ring structures, thus facilitating the phenomenon of appear- 
ance of spiral arms (Romero-Gomez et al. 2007), this parameter 
is not sufficient in order to characterize this phenomenon. The 
azimuthal displacement of the Lagrangian points is the most im- 
portant parameter. This result probably provides a dynamical ba- 
sis for understanding the reported failure of pure bar models to 
reproduce the inner spiral arms emanating at the ends of bars in 
both particle and hydrodynamical simulations of barred galaxies 
(e.g. Lindblad et al. 1996; Aguerri et al. 2001). 

The paper is organized as follows: Section 2 gives the form 
of the invariant manifolds in the models of Kaufmann and 
Contopoulos. We examine the manifolds a) when the spiral per- 
turbation is turned-on, and b) in an 'aligned model' version in 
which the whole non-axisyimnetric perturbation is artificially 
aligned to the bar. In case (a) the manifolds yield a spiral re- 
sponse, while in case (b) they yield a ring-like response. Since in 
all the above models the spiral perturbation is strong, we also ex- 
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amine two models corresponding to a 'mean' and 'weak' spiral 
amplitude, created by suitably varying the parameters of some 
of the original models of Kaufmann and Contopoulos (1996). 
We finally provide a theoretical justification of the importance 
of the azimuthal displacement of the Lagrangian points in the 
form of the invariant manifolds. Section 3 discusses the connec- 
tion between the theory of the invariant manifolds and the 'spiral 
response' models constructed via iterative methods. In particu- 
lar, we propose a method of constructing response models on the 
basis of populating by matter the manifolds generated by a 'pure 
bar' model. We also calculate response models via the method 
proposed by Patsis (2006) and discuss a number of morpholog- 
ical features of these models which find a straightforward ex- 
planation by the invariant manifolds. Section 4 summarizes our 
conclusions. 

2. Model and invariant manifolds 

2.1. Model 

The model of Kaufmann & Contopoulos (1996) consists of a 
number of potential/density terms representing various compo- 
nents of a baiTed-spiral galaxy. In particular we have: 
- A halo density term given by a Plummer sphere 



Phir) 



3M,, 



1 + 



2 ^ 

r 



-5/2 



- A disk surface density given by an exponential law 
2rf(r) = 2() exp(-erfr) . 

- A FeiTers bar with major axis aligned with the y-axis 



Pb(x,y,z) = 



105Mi 



32nabc 



1 



y 



X 



,2^ 



a > b > c. 



(4) 



(5) 



(6) 



Table 1. Parameters of models A, B, C (from Kaufmann & 
Contopoulos 1996). The units are km^s^^kpc ' for A, kpc ' for 
fij, Bd, K\ and K2, kpc for r\, r^. A, a, b, c, and b/,, lO'^M© for 
M/„ Mh, kms^'kpc"' for Q^, and M^/pc^ for Eq. Model A' has 
the same parameters as model A, except for the spiral ampli- 
tude A - 1000 and the pitch angle io - -9°. Model B' has the 
same parameters as model B except for A = 2500, Mh = 0.09, 
io = -8°, Qp = 30.5 kms-'kpc^^ 
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with 

2 2 2 

(p- b^ 

- and a spiral perturbation in the potential 
Vs(r, 0) - A(r)r exp(-e5r) cos2<I) 
where 

In(r/fl) 



0) = 

and 

A(r) = 



tan / 



lA-Ar 



if r > a 
if r < fl 



(7) 



(8) 



(l+tanh[^i(r-ri)])(l+tanh[^2('-2-r)])+A,.(9) 



The latter formula allows for a nearly constant amplitude A(r) 
between the inner and outer cut-off radii r\ and r2, while the am- 
plitude falls to a small value A^ beyond r2 or below r\ . We have 
in fact slightly modified Eq.© so that the potential becomes a 
smooth function of 9 r — a. Namely, we substitute ([8]) by the 
expression: 



1 + tanh(2(r - a)) 



ln(r/fl) 



tan I 



(10) 



The latter expression introduces a smoothing of the potential 
yielding a difference with respect to Eq.® which is 10% at the 
distance r - a + 0.5Kpc and only 2% at r = a + IKpc. 



A model is specified by a set of values for the parameters 
M/„ bk, Y.Q, e^. Ml,, a,b,c, 6,, A, A,., a^i, rj, K2, and r2, as well 
as the value of the pattern angular speed D.p. In Kaufmann & 
Contopoulos (1996), the parameters were adjusted so as to pro- 
duce three different self-consistent models which present some 
features of three real barred galaxies. The criterion for self- 
consistency was that the 'response density', i.e., the density ob- 
tained by the superposition of many orbits in the fixed potential 
should match as closely as possible the imposed density rep- 
resented by the above equations. The matching refers to a) the 
amplitudes of the surface density map on the disk plane, and b) 
the phases of the maxima of the bar and of the spiral arms, in the 
imposed and in the response models. 

The parameters for the three models are given in Table 1 . In 
the sequel we refer to these as model A, B, and C. The value 
of Ar for all three models, as well as the values of r2 and K2 for 
model C are missing from Kaufmann & Contopoulos (1996), 
where, however, it is noted that any (small) value of Ar, or of 
ri , r2 and ki , K2 does not influence the self-consistency (provided 
that r2 is at the end of the spiral arms). For consistency with the 
remaining models, we have set K2 - ki and r2 = 10 Kpc in the 
case of model C, and A,- = in all three models. 

Models A, B, C present some features of the galaxies 
NGC3992, NGC1073, and NGC1398 respectively. As discussed 
below, the bar-spiral strengths induced by the parameters of 
model C are quite untypical of barred-spiral galaxies, although 
still in the range allowed by observations. Thus, while the theory 
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Fig. 2. The ratio of non-axisymmetric forces, due to the bar or to the spiral arms, versus total axisymmetric force, as a function 
of the distance R along the x-axis (solid) or y-axis (dashed) for models A, B, and C (panels a, b, and c respectively). The vertical 
dashed lines mark the distance of the Li or L2 points in each case. 



of the invariant manifolds worked well in all three models, we 
discuss in detail models A and B, and only some exceptional fea- 
tures of model C interesting for dynamics. It should be pointed 
out that, while the choice of model parameters was partly based 
on observations (see Sects. 2, 3 of Kaufmann & Contopoulos 
1996), the so-obtained models are only rough representations 
of the referenced galaxies. For example, images of the galaxy 
NGC3992 (e.g. in the I band, TuUy et. al. 1996) indicate the 
presence of at least one more arm of amplitude comparable to 
the main bi-symmetric pattern. Images of the galaxy NGC1398, 
(e.g. in the R-band, Hammed & Devereux 1999) reveal the exis- 
tence of an inner ring structure which has no clear-cut separation 
from the main spiral structure. Such morphological features are 
not captured by the potential/density model given by Eqs. (|4]i- 
( fTOl ). Finally, the use of a n = 2 Ferrers bar model implies a 
steep drop of the bar force beyond the bar's limit which would 
be smoother in a n = or n = 1 model, and it does also not 
account for a rectangular-like outline that is observed in many 
real bars. 

These facts notwithstanding, the choice of potential param- 
eters and pattern speeds as in Table 1 ensures the existence of 
a self-consistent solution for the response density, a fact which 
would by no means be implied in an arbitrary choice of potential 
model. Although we do not make explicit use of the library of 
orbits of the final solution in the present paper, and also no guar- 
antee for the stability of the models is provided in Kaufmann 
and Contopoulos (1996), the self-consistency property suggests 
that the spiral arms found in these galaxy models can be stellar 
dynamically supported. This conclusion is independent of the 
theory of the invariant manifolds, thus the latter theory can be 
tested against this conclusion. 

The relative importance of the various non-axisymmetric 
components of the force with respect to the axisymmetric force 
vary with the distance from the center, as can be inferred from 
Figure 2. The bar contributes to the forcing by both an axisym- 
metric and a non-axisymmetric component. The axisymmetric 
component is found as the azimuthally averaged radial bar force 

FtarAr)^^f Fha,;Ar,9)de , (11) 
27r Jo 

where F),ar,r - (F bar,xX + F har^yj) I r . The total axisymmetric force 
is Fax{r) = Fdiskir) + Fhaio{r) + FharAr)- The non-axisymmetric 
bar force at a position r, 6 is the difference ^„on-ax - Ffenr - 



Fbar.r^r whcrc c, dcnotcs the unit vector in the radial direction. 
On the other hand the spiral force corresponding to the poten- 
tial can be all considered as non-axisymmetric since the spiral 
potential only has a cos 2<1) dependent term. Figure 2 shows the 
absolute ratios Fbar,non-axlFax and F.pirai/Fax as a function of the 
radial distance r, along two directions, i.e., along the bar's ma- 
jor (dashed curves) and minor (solid curves) axes, for all three 
models. 

In model A (Fig. 2a) the bar yields the dominant non- 
axisymmetric perturbation at all distances up to a zone around 
corotation (shown as a vertical dashed line at r = 6.11Kpc = 
Til ). The maximum amplitude of the non-axisymmetric bar force 
is 0.32, corresponding to a peak of the Ftarj curve atR ^ 4Kpc 
(in all the panels of Fig. 2 the innermost local maxima or min- 
ima of the curves Fhar,x, Fi,ar,\ al R < IKpc are artificial, due to 
the weakening of the axisymmetric forces which, for a Plummer 
sphere, are exactly equal to zero at R - 0). The inner width of 
the zone is found by the point where Fhar,non-ax = F spiral, which 
is at a distance r ^ 5.75Kpc. Beyond that distance, the spiral 
term dominates over the bar term, reaching a maximum ampli- 
tude equal to 0.21 with respect to the axisymmetric background. 
The oscillations of the spiral force beyond corotation are due 
to the logarithmic dependence of the argument <1> in (|7]l on r, a 
fact causing successive maxima and minima of the spiral force 
at successive periods of length In of the argument The first 
maximum, around corotation, is the most important. The width 
of the oscillation from this maximum to the next defines an ap- 
proximate value of the radial wavelength of the spiral density 
wave, which is Ar ^ IKpc. 

In model B (Fig. 2b) the maximum amplitude of the non- 
axisymmetric bar force reaches the value 0.75 (for Fbar,v at 
R ^ 2Kpc), implying that the bar is quite strong inside corota- 
tion. Nevertheless, even in this galaxy the spiral force becomes 
dominant over the bar's non-axisymmetric perturbation around 
and beyond corotation. The zone around the first maximum of 
the spiral force defines a radial wavelength of the spiral den- 
sity wave Ar ^ IKpc. The first peak of the spiral force is again 
found to be at a distance very close to the corotation radius 
r = ru - 3.43 and the amplitude of this peak is 0.35. 

Finally, in model C (Fig. 2c) the spiral perturbation near and 
beyond corotation reaches such a high amplitude (maximum = 
0.57), that it becomes even stronger than the maximum ampli- 
tude of the bar's perturbation (^ 0.3 for Fhar,y at R ^ 3.5Kpc) 
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Fig. 3. The Qs vs Qb values in observed SB galaxies (gray points 
from the Buta et al. (2005) sample) and our models A, A', B, B', 
C (black dots). 



which takes place well inside corotation. Furthermore, the spi- 
ral arms are tightly wound (the radial wavelength is estimated 
as Ar ^ O.SKpc), and the spiral arms extend to cover about one 
azimuthal period 2n. Thus, model C is exceptional and will not 
be discussed in detail in the sequel. Only a feature of this model 
interesting for dynamics is discussed in Sect.2.3. 

The bar-spiral amplitudes of models A and B define strongly 
nonlinear models, which are above the average but inside 
the range of bar-spiral strengths found by recent observations 
(Laurikainen & Salo 2002; Buta et al. 2005). In the latter works 
the maximum value of the ratio of the tangential force versus 
the radial force for the bar and spiral components is denoted 
by Qh and Qs respectively. The average values for SB galax- 
ies in the sample of Buta et al. (2005) are < Qb >- 0.29 and 
< Qs 0.17. Model A yields Qb = 0.3 and Qs = 0.21. 
Model B yields Qb - 0.75 and Qs - 0.35, both values being 
a factor 1.3 larger than the specific estimates reported by Buta 
et al. (2005) for the galaxy NGC1073 {Qb = 0.56, g, = 0.26), 
to which model B is associated. In order to have a more repre- 
sentative sample of models in which the theory of the invariant 
manifolds is to be tested, two 'weak' models. A' and B', are also 
considered, which were created by varying some parameters of 
models A and B. In model A' the spiral amplitude is A = 1000, 
i.e. half the value of model A (see table I). In model B' the bar 
mass is Mb - 0.09 and the spiral amplitude A - 2500. Since 
these changes are rather arbitrary, there is no guarantee of self- 
consistency of the new models. However, a rough criterion of 
self-consistency (subsection 2.2) can be established if the pitch 
angle is also slightly varied in both models (/q = -9" in model 
A' and /q = -8° in model B'). The resulting Qb and Qs values are 
Qb = 0.3, Qs = 0.125 for model A' and Qb = 0.5, g, = 0.15 for 
model B'. The spiral strengths of models A',B' are well below 
the average of the Buta et al. sample for SB galaxies. Thus, the 
mean values < Qb > and < 2^ > of the four models A,B,A',B' 
become both representative of the average values found in the 
observations. On the other hand, the value Qs - 0.57 of model C 
is untypical although still in the range of the observations (Figure 
3). 



2.2. Phase portraits and invariant manifoids 

The first result of the analysis of the invariant manifolds can now 
be demonstrated with the help of Figs. 4 to 6. Figure 4a shows 
the phase portrait (surface of section (6*, pg) coiTesponding to the 
apocentric positions r - Q, pr < 0) in the case of model A, 
for a value of the Jacobi constant Ej = -1.91 x 10^, which is 
close to the value Eju = -1.915 x 10^. Figure 4c shows the 
same portrait for Ej = -1.91 xlO^inaso called ' aligned spiral' 
version of model A in which the angle O in Eq.dTOt is replaced 
by 6 throughout the whole radial extent of the spiral arms. This 
means to artificially 'align' the spiral arms as extensions of the 
bar along the latter's major axis. By this way we measure the 
effect of only increasing the amplitude of the non-axisymmetric 
perturbation on the form of the invariant manifolds, while in the 
original model the manifolds are aff'ected both by the strength 
of the non-axisymmetric perturbation and by the azimuthal dis- 
placement of the unstable Lagrangian points with respect to the 
bar's major axis. 

In Figs.4a,c the points marked PLl, PL2 correspond to 
the fixed points of the PLl or PL2 short-period orbits which 
are close to the positions of the unstable equilibria Li, L2. 
Furthermore, the thick dots show the intersection of the unstable 
manifolds "Wp^^ and 'Wp^2 with the surface of section. In or- 
der to facilitate the reading of these diagrams, we note that, for 
Pr = (apsides), beyond some radius ro ^ IKpc, Eq.([TJ yields 
that r increases nearly monotonically with pe in all azimuthal 
directions of a fixed angle 6 (a small reversal of this monotonic 
relation, due to the non-axisymmetric potential terms, is only 
observed at angles 6 ^ jt/2 and in a small interval of radii, of 
width Ar < O.lKpc around r = 5Kpc; the monotonic relation is 
re-established after crossing this interval). Thus, in Figs.4a,c the 
semi-plane of the phase portrait with pg > pgxi means apocen- 
tric positions outside corotation, while pg < pe n means apoc- 
entric positions inside corotation. Note that in this and in all sub- 
sequent plots of phase portraits the values of pg are normalized 
with respect to the value a^Qp, corresponding to the angular mo- 
mentum in the rest frame of a circular orbit at a radius r - a. 

The main remarks about the comparison of the two phase 
portraits are now the following: 

- In both portraits chaos is pronounced inside corotation (for 
pg < pg,u), and the domain of inner invariant KAM curves is 
deeply inside the bar (at values of pg about or below 0.25). Such 
extended chaotic domains are responsible for the termination of 
the bar. 

- Outside corotation (for pg > pg^i), a layer of outer KAM 
curves has been destroyed in both portraits. This is caused 
mainly by the growth of the chaotic layer around the unstable 
-6/1 periodic orbit, which produces a resonance overlap with the 
chaotic layer of the PL 1,2 unstable periodic orbit (negative signs 
indicate a resonance outside corotation, for which the motion is 
retrograde in the azimuthal direction). As a result, the chaotic 
domain extends up to values of pg ^ 1.5, and the first rota- 
tional KAM curves appear a little inside the -4/1 resonance. The 
domain around the outer Lindblad resonance is almost entirely 
filled either by rotational KAM curves or by 'resonant' curves 
around the -2:1 stable periodic orbits. The islands of stability 
of the -2:1 resonance have a larger width in the 'aligned spiral' 
model (Fig. 4c) because by aligning the spiral perturbation the 
amplitude of the total non-axisymmetric perturbation increases 
effectively at large distances from corotation (the width of reso- 
nances scales as a power-law of the non-axisymmetric perturba- 
tion (see e.g. Contopoulos 2002). 



p. Tsoutsis et al: Invariant manifolds and spiral arms 



7 




Fig. 4. (a) Phase portrait near corotation (Ej - -1.911 x 10^) in the case of model A. The thick gray points mark the position of 
the fixed points of the PLl and PL2 orbits. The thick dark lines are the invariant manifolds "Wp^^ and "W^^^- (t') Projection of the 
manifolds of (a) in the configuration space. The gray spiral curves correspond to the density maxima of the imposed spiral arms of 
the full model. The curves marked LC4,LC5 and KAM are commented in the text, (c), (d) Same as (a), (b) in the 'aligned spiral' 
{<S> - 0) version of model A, in which the minima of the spiral potential term are aligned with the bar, and Ej - -1.91 x 10^. 



- The white circular domains devoid of points, embedded in 
the chaotic sea of both portraits, correspond to prohibited do- 
mains of motion, for r = and for the selected values of the 
Jacobi constant. Such domains exist when Ej < E (equal to 

= -1.873 X 10^). 

- The most important difference between the two portraits 
is that in the case of the true spiral term turned on (Fig. 4a) the 
prohibited domains lose their azimuthal symmetry with respect 
to the values 6 = (equal to 2n{mod2n)), or = n referring 
to the positions of the stable Lagrangian points L4 and L5. Such 
a symmetry is perfect in the aligned spiral case (Fig. 4c). The 
limiting boundaries of the prohibited domains are denoted by 
'LC4', 'LC5' in Fig.4. The azimuthal deformation of the pro- 
hibited domains corresponds to an azimuthal deformation of the 
associated banana-like prohibited domains appearing in the con- 
figuration space (i.e. the disk plane). These domains are simi- 
lar but should not be confused with the domains delimited by 
the zero velocity curves of the effective potential in the rotating 
frame, i.e. V,ff = V(r, 6) - Q^rV2. 



It turns out that the azimuthal deformation of the limiting 
boundaries LC4 and LC5 is a crucial diff'erence related to the 
morphology of the spiral arms. In both Figs.4a,c the manifolds 
are well developed inside and outside corotation (below and 
above PLl). In particular, the invariant manifolds form conspic- 
uous lobes and foldings which are typical of systems having a 
large degree of the so-called homocUnic chaos. However, when 
the manifolds are plotted in the configuration space (Figs.4b,d), 
the azimuthal deformation of the limiting boundaries LC4 and 
LC5 in the truly spiral model (Fig. 4b) causes the manifolds to 
be also azimuthally deformed in a way so as to closely support 
the imposed spiral arms up to an azimuth 9' ^ 37t/4 measured 
clockwise from either Li or L2S. The end of this agreement is in 
two regions (marked A and A') in which the manifolds "Wp^^ 
and 'Wp^r, lose contact from the locus of maxima predicted by 
Eq.®. Beyond this distance, the manifolds exhibit a typical be- 
havior called a 'bridge' in our previous works (Voglis et al. 2006, 
Tsoutsis et al. 2008). Namely, the manifold 'Wp^^ forms an in- 
ner spur consisting of a number of lobes connecting segments 
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Fig. 5. Same as in Fig.4, but for model B. The Jacobi constant is Ej = -2.98 x 10'' in (a) and (b), and Ej = -2.97 x W in (c), (d). 



of the manifold from A to their continuation, which starts in the 
neighborhood of the other periodic orbit, PL2. From that point 
on the oscillations of the manifold 'Wp^^ start supporting both 
the spiral arm emanating from L2 and the chaotic layer marking 
the border of the bar. When the manifold is calculated for even 
longer lengths, we find that the higher order lobes of the mani- 
fold from Li, after reaching the neighborhood of L2, come again 
close to Li in directions nearly parallel to those formed by the 
low order lobes of the same manifold. This causes an enhance- 
ment of the density close to L\ along the direction of the imposed 
spiral arms, i.e., the manifolds support the self-consistency of the 
spiral arms. 

On the other hand, the manifolds of the 'aligned spiral' ver- 
sion of model A (Fig.4d), calculated up to a length comparable to 
that of the manifolds of Fig.4b, show no support of a spiral struc- 
ture, but only yield a thick ring-like structure. The thickness of 
the manifolds of Figs.4b,d is determined by the degree of chaos 
in Figs.4a,c. The degree of chaos is determined by the amplitude 
of the non-axisymmetric perturbation. This is expected from dy- 
namical systems theory, since the overlapping of resonances, 
which is the main source of production of chaos, depends on 
the width of the different resonant layers near corotation, which, 
in turn, depends on only the amplitude of the perturbation. 



A more elaborate analysis (Sect. 2.3) shows that, while the 
outermost radial limit of the invariant manifolds is posed by the 
existence of absolute barriers, i.e. rotational KAM tori (marked 
'KAM' in Fig.4b), more stringent limits are practically posed by 
partial barriers, i.e. cantori, which hmit the diffusion within a 
chaotic zone. Provided these limits, the azimuthal deformation 
of the invariant manifolds is the crucial factor for the production 
by them of response spiral arms. This, in turn, is determined by 
the form of the hmiting boundaries LC4 and LC5. The theoreti- 
cal derivation of these boundaries is given in (Sect.2.3). 

Figure 5 shows the same phenomena in the case of model B. 
The qualitative resemblance between Figs.4a,b and 5a,b is obvi- 
ous, although the azimuthal deformation of the limiting bound- 
aries LC4 and LC5 is more pronounced in Fig. 5b than in Fig.4b. 
Also in this model the manifold exhibits a bridge starting at an 
angle ^ 3;r/4 clockwise from Li or L2 (points A, A'), as well 
as inner spurs connecting segments of it along both spiral arms 
and along the border of the bar. Another feature of Fig. 5b is that 
the inner branch of the invariant manifold (inside the bar) is de- 
veloped in a domain occupying about one fourth of the total ex- 
tent of the bar. This imphes that a substantial part of the bar in 
the domain near corotation is supported by chaotic orbits. In fact 
the non-axisymmetric forcing in model B is much stronger in- 
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Fig. 6. (a,b) Same as in Fig.4a,b but for the model A', and Ej = 
Ej = -2.77 X lO'*. 

side corotation than in model A, a fact causing the destruction 
of all the inner KAM curves down to pg = (Fig. 5a). Such 
a type of chaos may lead to a number of observational conse- 
quences, photometric and kinematic, a list of which have been 
enumerated by Grosb0l (2003). Finally, the azimuthal deforma- 
tion of the maxima of the spiral term with respect to the bar's 
major axis also turn out to be the crucial factor for the produc- 
tion by the manifolds of response spiral arms. In fact, by com- 
paring Figs.5a,b with the respective figures in the 'aligned spi- 
ral' version of model B ( Figs.5c,d) we see that the manifolds in 
the latter case present some asymmetry as well as a large thick- 
ness, due to the high value of the non-axisymmetric perturbation, 
but they still largely deviate from the spiral pattern (gray locus), 
which was closely followed by the manifolds of the non-aUgned 
model (Fig.Sb). 

Figure 6 shows the phase portrait structure near corotation in 
the 'weak spiral' models A' (Fig. 6a) and B' (Fig.6b). The invari- 
ant manifolds 'Wp^j and ^p^j plotted, and the coun- 
terparts of these plots in the configuration space are shown in 
Figs. 6b and 6d respectively. As expected, in both models chaos 
is considerably reduced with respect to the strongly nonlinear 
models A,B, and it is only limited in a narrow zone in the coro- 
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-1.91 X 10^. (c,d) Same as in Fig.5a,b but for the model B' and 



tation region. Further away, the phase space is filled by invariant 
tori which occupy most of the phase space volume already at the 

-4; 1 resonance. 

The Qs value of model A' is Qs = 0.125, and this is its only 
diflFerence with respect to model A, which has Qs = 0.21. The 
thickness of the invariant manifolds is thus reduced with respect 
to the thickness of the manifolds of model A (compare Figs. 6b 
and 4b). However, the azimuthal deformation of the manifolds 
is still large enough to fit the locus of the imposed spiral arms 
up to an angle = n/2 clockwise from L\ or L2, i.e. the man- 
ifolds support quarter turn spiral arms. In the case of model B' 
(Fig.6d) we have Qs = 0.15, which is close but still below the 
average value < g j >= 0.17 of the Buta et al. (2005) sample. At 
the value Qs = 0.15 the region of homocUnic chaos formed by 
the lobes of the manifolds near Li or L2 is already well devel- 
oped, and the 'inner spurs' are clearly distinguishable. In fact, in 
model B' a small adjustment of the pattern speed (Q.p - 30.5 in- 
stead of 32.5 Km/sec/Kpc, as was in model B) yielded the best fit 
of the invariant manifolds to the imposed spiral arms. Such a fit 
can now be considered as a rough criterion of self-consistency. 
We see that the formation of bridges in the manifolds of Fig.6d 
result in that the higher order lobes of the manifold make oscil- 
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lations which enhance the density along the manifolds' unstable 
directions all the way from Li or L2. Thus the manifolds sup- 
port, again, the imposed spiral arms in a self-consistent way. In 
fact, since the overall thickness of the domains covered by the 
invariant manifolds increases, the manifolds support the spiral 
structure up to an angle larger than n/2, i.e. along a length 
larger than in model A'. 

In conclusion, there is a clear morphological continuity of 
the structures produced by the invariant manifolds, from rings 
to quarter turn spirals, and then to fully developed spiral arms, 
as the value of Qs increases. The azimuthal displacement of the 
unstable Lagrangian points is responsible for the manifolds pro- 
ducing a spiral-like response, while the increase of the amplitude 
Qs extends the support of the spiral structure to higher angles 
0' . Such a morphological continuity of the invariant manifolds 
is suggestive of it being a real morphological feature of barred 
galaxies. 

2.3. Theoretical modeling 

A first estimate of the azimuthal deformation of the limiting 
boundaries LC4 and LC5 can be done by calculating the az- 
imuthal displacement of the unstable Lagrangian points Li 2 
when Vspimi is turned on (Figs.7a,b, schematic). This can be 
judged from the form of the equipotential curves (called here- 
after the 'curves of zero velocity', CZV) in the rotating frame, 
i.e., the level curves of 



Veff^V(r,0)-Qpr^/2^Ej 



(12) 



For a logarithmic spiral, the azimuthal shift of Li is given by 



A9 = 



1 



tanio 



In^ 



(13) 



where AO is taken clockwise from the bar's major axis. The az- 
imuthal deformation of the CZVs takes place in a narrow zone 
around corotation (Fig. 7b), while, for r small, the curves remain 
practically aligned to the bar This is because the spiral potential 
term ^ (with <t) given by Eq.lfTOll) adds practically a cos 26 con- 
tribution to the bar's term at radial distances r < a. The stable 
Lagrangian points L4 5 are shifted counterclockwise with respect 
to the bar's minor axis, although by a smaller angle than that of 
the unstable points L\, L2. 

The equipotential curves of Vefj do not provide the strictest 
limit of allowed motions on the apocentric surface of section 
(0, pe), or (6, r), for pr - 0, pr < 0. The form of the limit- 
ing curves in the configuration space is shown schematically in 
Fig. 7c (curves LC4 and LC5, encircling the stable Lagrangian 
points L4, L5). These curves are derived by the request that, for 
any fixed angle and angular momentum pg, the curve of the 
function 



pI 

V(r,0,pe)^-^-Q„Pe + V(r,0) 



(14) 



be tangent to the line V - Ej at the point where V - Ej has 
a global minimum. The minimum is calculated by the root r - 
rm(d,Pe) of 



dV 



r' dr 



= 



(15) 



The hmiting values of pg, for fixed 0, are then found via the roots 
for pg of 

V(r,Mp8),0,pg)^Ej . 



The limiting curves LC4 and LC5 in the configuration space, 
derived from Eqs.(fT4li and ( fTSl ) are outside the limiting curves 
provided by the CZVs defined through Eq.(fT2b. This is due to 
the fact that these limits now refer to r = 0, while in the case of 
Eq.dTSb the total velocity in the rotating frame is equal to zero 
(details are given in Appendix A). By virtue of these facts, the 
curves LC4 and LC5 exhibit also a spiral-like azimuthal defor- 
mation, which is even more pronounced than that of the CZVs. 

The gray domain in Fig. 7c shows the permissible apocentric 
positions of the orbits under a fixed value of the Jacobi constant 
close to the corotation value. The interior gray domain between 
LC4 and LC5 roughly marks the extent of the bar. On the other 
hand, the positions of the PLl and PL2 points are in very narrow 
strips of permissible apocentric positions separating the inner 
right part of the LC4 curve from the outer right part of the LC5 
curve and vice versa. The invariant manifolds Wp^j j emanating 
from these points necessarily follow the narrow strips leading to 
the outer gray domain, thus they yield locally the form of spiral 
arms. 

The evolution of the invariant manifolds further away from 
PLl or PL2 is determined by the resonant structure in the outer 
corotation zone. The existence of many resonances accumulat- 
ing in a narrow range of distances near the corotation radius 
causes a chaotic layer in this region, formed by the mechanism 
of resonance overlap. Figure 8 makes a zoom to the phase por- 
traits of the models considered, in order to demonstrate the rel- 
evant phenomena. Figs.8a,b are zooms to the phase portraits of 
Figs.4a, and 5a, referring to the models A and B. In both cases 
we find a chaotic layer extending up to pg ^ 1.5, which is de- 
limited by a rotational KAM torus at pg = 1.5 (marked KAM). 
This torus is just below the -4/1 resonance, which is stable at 
the value of the Jacobi constant Ej - -1.911 x 10^. On the 
other hand, most islands of stability of resonances -m/l, with 
m > 4, which are closer to corotation, have been destroyed. 
Their destruction is followed by the destruction of KAM tori 
with irrational rotation numbers, which, according to the stan- 
dard theory, are transformed into cantori. Such cantori limit the 
chaotic flux through their gaps, and this fact causes some stick- 
iness in a zone very close to the PLl and PL2 fixed points. 
Stickiness phenomena of this type have been explicitly demon- 
strated and studied in simple models of the dynamical systems 
theory (see e.g. Efthymiopoulos et al. 1997; Contopoulos et al. 
1999; Contopoulos & Harsoula 2008), and they have also been 
observed in our A^-Body simulations of barred-spiral galaxies 
(Tsoutsis et al. 2008). The main outcome of these studies is that 
cantori in a large chaotic sea act as partial barriers slowing down 
considerably the escape (or diffusion) of the chaotic orbits with 
initial conditions along or near an invariant manifold. 

In the case of the manifolds "W^^^ , plotted in Figs.8a,b, 
which are calculated from 21 iterations of an initial segment 
of length ds - 10 close to PLl or PL2, we see that the 
manifolds fill only partially the chaotic domain up to the torus 
marked KAM. The inner dark region covered by the first iter- 
ations of the invariant manifolds defines a domain called 'in- 
ner stickiness zone', the projection of which in the configura- 
tion space is shown schematically as a dark gray domain in 
Fig. 7c. The remaining part of the chaotic domain up to the curve 
marked 'KAM' corresponds essentially to the light gray domain 
of Fig. 7c. This domain is eventually covered by the invariant 
manifolds after a very large number of iterations. For example, 
the manifolds of Figs.8a,b have not yet reached the curve KAM 
after about 50 iterations, which in both models correspond to 
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Fig. 7. The curves of zero velocity (equipotential curves of the effective potential Ve//) shown schematically when the non- 
axisymmetric perturbation is (a) aligned to the bar, and (b) has a spiral component, (c) A schematic figure of the basic theory 
according to which the structure of the chaotic phase space in the corotation region supports of trailing spiral arms. The spiral arms 
are confined in an inner stickiness zone (dark) delimited by cantori which are the remnants of rotational KAM tori outside corota- 
tion. Inside the bar, the same chaotic zone is delimited by rotational KAM tori (regular bar domain). The chaotic zone between the 
curves LC4 and LC5 and the boundary of the stickiness zone has a spiral shape, due to the deformation of the limiting curves LC4 
and LC5 caused by the azimuthal displacement of LI and L2 with respect to the bar's major axis. 



about 60 pattern rotation periods. Thus, during all this time in- 
terval the manifolds support a spiral structure. 

The stickiness phenomena keep playing a significant role 
even when the spiral perturbation is pushed to untypically high 
values. For example, in model C, (Figs.8c,d) the spiral strength 
is Qs - 0.57, and under such a high value all the rotational KAM 
curves are destroyed, at least up to the outer Lindblad resonance. 
Then, while in principle there is no absolute barrier to chaotic 
diffusion up to very large distances from corotation, a plot of the 
invariant manifolds (Fig. 8c) shows that these manifolds exhibit 
again stickiness phenomena, and they practically remain con- 
fined for very large times below the -4/1 resonance (which is 
still stable at the value of the Jacobi constant /sy - -2.9Y1 x 10^, 
yielding four tiny islands embedded in the large chaotic sea of 
Fig. 8c near the level pe - 1.6. This results in that the mani- 
folds in the configuration space (Fig.8d) yield the form of tightly 
wound spiral arms. In this particular example, the domain cov- 
ered by the invariant manifolds practically coincides with the 
'inner stickiness domain' of Fig. 7c. 

3. Spiral arms as the response of invariant 
manifolds to bars 

One immediate consequence of the analysis of the previous sec- 
tions is that one cannot induce the morphology of the spiral arms, 
corresponding to a particular morphological type of bar, by cal- 
culating the invariant manifolds of the PLl and PL2 families in 
only a pure bar potential. In fact Yspirai is most important near 
corotation and it must be taken into account self-consistently 
in all studies related to the morphology of the spiral arms via 
the calculation of invariant manifolds. This result is in agree- 
ment and probably provides a dynamical basis for understand- 
ing the results of both particle and hydrodynamical simulations 
(Lindblad et al. 1996; Aguerri et al. 2001) which have reported 
the inefficiency of simulations of pure bars to reproduce a spiral 
structure. 

On the other hand, the theory of the invariant manifolds sug- 
gests that the spiral arms are linked dynamically to the bar A 
plausible scenario for establishing such a link is one in which 



the bar initiates the process of a spiral response, which is then 
enhanced self-consistently by the growing contribution of the 
spiral potential. 

The question addressed below is whether, on the sole knowl- 
edge of the gravitational potential and pattern speed of a pure 
bar model, such a process can be modeled via the theory of the 
invariant manifolds. We construct such an iterative 'spiral re- 
sponse' model as follows: 

i) We first calculate the invariant manifolds produced by the 
pure bar model. 

ii) We assume that the invariant manifolds produced by the 
pure bar alone 'trigger' the formation of a ring-like or spiral pat- 
tern by attracting matter along the invariant manifold. The den- 
sity gradient along the manifold cannot be uniform, since (a) the 
speed of chaotic diffusion is smaller close to L\ or than far 
from these points, and (b) the higher order lobes of the invariant 
manifolds return to the neighborhood of L\ and L2 (Sect.2.2). 
In order to model the mass distribution along the invariant man- 
ifolds, we consider a number A^^ of small Plummer spheres, of 
mass nip and softening radius r^, placed along the invariant man- 
ifold of the pure bar case (Fig. 9a) with a linearly decreasing mass 
from L\ or L2 counterclockwise, namely the mass of the i-th par- 
ticle is given by: 



m„{i) 



2M, 



2M, 



(16) 



with i — 1,2, ...,Npl2 from L\ to L2 counterclockwise and simi- 
larly along the symmetric manifold from L2 to Li . Ms is an esti- 
mate of the total mass on the spiral arms (see below). The choice 
of a linear mass decrease as in (fTST i along the response spirals is 
rather arbitrary and it does not follow directly from the theory 
of the invariant manifolds. However, it does capture the essen- 
tial feature that the density of points should in general decrease 
along the unstable manifold as we recede from the unstable pe- 
riodic orbit. 

In the simulation of Fig. 9, we start from the invariant man- 
ifolds of the pure bar version of model A, and set rp - 0.5, 
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Fig. 8. The resonant phase space structure in the corotation region in the cases of (a) model A with Ej - -1.911 xlO^, (b) model B 
with - -2.98 x 10"*, and (c) model C with isy = -2.912 x 10^. The invariant manifolds of the PLl and PL2 orbits are over-plotted 
(thick lines). The bold rotational KAM curves near pg - 1.5 in (a) and pg - 1.57 in (b) are very close to the inner boundary of 
the regular domain beyond corotation. In (c) there is no such boundary at least up to the outer Lindblad resonance. The projection 
of the rotational curves of (a) and (b) marked KAM, as well as of the curves LC4 and LC5, in the configuration space, yields the 
respective outer closed curves in Figs. 4b, 5b. (d) Same as in Fig.4b but for the model C, and Jacobi constant Ej - -2.912 x 10^ 



Np - 30, and the mass of each particle fixed so that the total 
mass of all the particles is equal to M, given by 

'»r2 r'ln 



^ Jr, Jo J-o, 



\ps(r,8,z)\rdrd6dz 



(17) 



where the quantity ps(r, 0, z) is an approximate expression for the 
density perturbation corresponding to the potential (|7]l given by 
the WKB ansatz: 



Ps(r,0,z) 



mm\ + 1) 



V,(r, 0)sech2+l*'^l(z/A) 



(18) 



where k - 2/(rtani), according to the formula given by 
Vandervoort (see Contopoulos and Grosb0l 1988; the used val- 
ues of A are given in Table 1). Thus, Ms represents an estimate 
of the total mass contained in the spiral arms of model A. A gray 
scale plot of the total surface density corresponding to the above 
mass distribution is shown in Fig. 9b. 



iii) The gravitational potential is calculated anew taking 
into account the 'response' potential produced by the Plummer 
masses positioned as described in step (ii). This also yields a new 
position of the Lagrangian points Li, L2 as well as a new form 
of their invariant manifolds. The procedure should be repeated 
until convergence of the positions of the Lagrangian points and 
of their manifolds towards a final form is obtained. In practice, 
we find that just one iteration suffices to obtain invariant mani- 
folds which fit the imposed spiral perturbation up to an azimuth 
(ff ^ 7r/2). In particular, the manifolds of Fig. 9c cover a much 
thicker radial domain than those of the initial 'pure bar' model 
(Fig. 9a), and they also exhibit an azimuthal deformation follow- 
ing from the non-uniform distribution of mass shown in Fig. 9b. 

The present method can in principle be used to explore the 
sequence of different morphological types of spiral arms ob- 
tained by the method of the invariant manifolds applied to a 
family of pure bar potentials, when the bar parameters, i.e., bar 
strength and pattern speed are altered. Such a type of study 



p. Tsoutsis et al: Invariant manifolds and spiral arms 13 




Fig. 9. A 'spiral response' model based an initially 'pure bar' version of model A. (a) Manifolds of the pure bar case for Ej - 
— 1.90x 10^ (very close to the Ejl\ value. These manifolds are populated by 'point masses' represented as Plummer spheres centered 
at the positions indicated by the the gray thick dots (see text for details), (b) Gray scale mapping of the surface density produced by 
the mass distribution of (a), (c) Response manifolds calculated in the new potential corresponding to the mass distribution of (a). 




Fig. 10. Unstable invariant manifolds of the PLl and PL2 orbits in a pure bar version of model A in which the bar's mass and pattern 
speed are altered with respect to the reference values Mj, = 1.5 x lO^^Mg, Q.p = 43.6Km/sec/Kpc so as to keep corotation at a fixed 
distance, (a) M; = 2Mb, ^'p = 1.260p, Ej = -2.19 x 10^ (b) M; = 3Mb, = 1.38Qp, Ej = -2.40 x 10^. 



was undertaken by Romero-Gomez et al. (2007), without tak- 
ing, however, into account the effects of the self-consistent spiral 
response. For our adopted models and parameters, we find that 
when this efi'ect is ignored, a kind of spiral response can only 
be produced by the low-order lobes of the invariant manifolds, 
but i) the bar strength has to be pushed to very high values, and 
ii) the spiral form disappears when the manifold is calculated 
for longer length. An example is given in Fig. 10, referring again 
to the 'pure bar' version of model A. The invariant manifolds 
are calculated after altering the bar parameters with respect to 
the values given in Kaufmann & Contopoulos (1996). Figure 10 
shows two different choices of parameters, in which the bar is 
given a mass equal to (a) two and (b) three times the mass value 
in Kaufmann & Contopoulos (1996). The pattern speed was al- 
tered accordingly so as to keep corotation at a fixed distance. 
As the bar amplitude increases, the first lobes of the invariant 
manifolds become more open and they yield a gradual transi- 
tion from a ring (Fig. 10a) to a spiral pattern (Fig. 10b). However, 



the spiral pattern in Fig. 10b does not fit the self-consistent spiral 
pattern of model A. Furthermore, the change of the bar strength 
required in order to produce this result constitutes a large devia- 
tion from the parameters of the original model which guaranteed 
self-consistency. This can be probably improved by choosing a 
lower n-value of the Ferrers bar and/or including the effects of 
a more rectangular bar outline. 

A different methodology to produce response models of 
barred-spiral galaxies has been proposed by Patsis (2006). As 
already emphasized in Tsoutsis et al. (2008), the theory of the in- 
variant manifolds explains many features of such response mod- 
els. Patsis (2006) considers an ensemble of particles with initial 
conditions on circular orbits of the axisymmetric part of the po- 
tential, placed uniformly on the disk at radii corresponding to 
Jacobi constants up to the L4 value of the full potential. Then, 
by softly introducing the non-axisymmetric part of the potential 
(the transition time is a few pattern periods), the particles' distri- 
bution changes due to both the adiabatic change of the potential 
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Fig. 11. Model A for the value of the Jacobi constant Ej - -1.911 x 10^. (a) The two gray horizontal lines define a zone of initial 
conditions around the circular orbit of the unperturbed (axisymmetric) potential. The thick black dots are the images of these initial 
conditions after one iteration of the Poincare map. The background thin points show the underlying phase space structure, (b) The 
first image of the Poincare map of the initial conditions specified in (a) as it appears in the configuration space. After only one 
iteration this set of points already acquires the form of the invariant manifolds, i.e., it yields trailing spiral arms. The projection of 
the tenth iterates of the same points on the surface of section and on the configuration space are shown in (c) and (d) respectively, 
(e) Same as in (d) when the initial conditions are taken in a zone four times as wide as the zone of (a), (f) A fuU 'response model' 
for the same galaxy produced according to the methodology suggested in Patsis (2006). The black curves represent the minima of 
the imposed spiral potential. 
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and phase mixing. As a result, the particles finally settle to orbits 
supporting both the bar and the spiral arms. 

Inside the bar, the particles of such response models are cap- 
tured mostly in stable resonances belonging to a branch or bifur- 
cation of the Xi family, which is the continuation of the family of 
circular orbits of the axisymmetric potential on which all the ini- 
tial conditions lie. However, near corotation and beyond chaos is 
prominent. In this case the theory of the invariant manifolds ex- 
plains the capture of the particles in orbits supporting the spiral 
arms. In our specific models, one can see this effect by taking 
initial conditions (6, pe) on the surfaces of section such as those 
of Figs. 4a and 5a such that these initial conditions a) belong to 
a chaotic domain of their respective surfaces of section, and b) 
determine circular orbits when all the non-axisymmetric part of 
the potential is turned off (the monopole contribution of the bar, 
obtained by averaging radial forces with respect to all possible 
azimuths, is taken into account in this calculation). The locus of 
all these initial conditions is a straight Une pg = const on the 
surface of section for a given Jacobi constant. In practice, we 
take a narrow zone of some width around such a line, i.e., allow 
also for a small value of the epicyclic action around the circular 
orbits of the axisymmetric model. 

Figure 1 1 shows the result of ruiming these initial conditions 
in the case of model A, with the non-axisymmetric part of the 
potential being turned on from the start, and for one iteration 
(Figs. 1 la,b, corresponding to 1 .2 pattern rotation periods), or ten 
iterations (Figs.llc,d, corresponding to 12 pattern rotation peri- 
ods), forthe value of the Jacobi constant = -1.91 1 x 10^. The 
initial conditions correspond to a zone in Fig. 11a between the 
two gray horizontal lines, that is we take 10000 points uniformly 
distributed within the intersection of the zone with the permissi- 
ble domain of motion. The first iterates of these points are shown 
with dark thick points in the same figure, while the large number 
of small dots illustrate the overall structure of the phase portrait 
at the chosen value of Ej. The corresponding figure in config- 
uration plane is shown in Fig. 1 lb. An obvious conclusion from 
Figs. 1 la,b is that, already after one iteration, the initially straight 
zone of initial conditions in the surface of section is deformed so 
as to closely follow a pattern induced essentially by the invari- 
ant manifolds "W^^j ^ (compare Figs. 1 la,b with Figs.4a,b). This 
phenomenon is repeated at subsequent iterations, so that after 
ten iterations (Figs. 1 1 c,d) the set of all points describes a pattern 
nearly coinciding with that of the manifolds IVp^^ o^^r ^ 
larger length of the latter. This picture does not change qualita- 
tively if the zone of initial conditions is taken to have a width 
four times as large as in Fig. 11a. The resulting response of the 
orbits in this case is shown in Fig. lie, and the association of this 
with the dynamics of the invariant manifolds is still quite clear. 

This behavior follows from a well known 'mixing property' 
of chaotic dynamical systems, namely that any small and com- 
pact ensemble of initial conditions embedded in a large chaotic 
domain is deformed in subsequent Poincare mappings so as to 
follow the unstable invariant manifolds of the main families of 
unstable periodic orbits located in the chaotic domain, while pre- 
serving its measure in the same domain (see Contopoulos & 
Harsoula (2008) and references therein for a detailed exploration 
of this phenomenon in the case of the 'standard map'). In the 
case of barred galaxies, this type of response of the chaotic or- 
bits to the invariant manifolds generates patterns such as those 
of Figs. 1 lb,d i.e., spiral arms. 

It should be stressed, however, that Figs. lib to lie were 
calculated with the full non-axisymmetric potential, i.e., Vhar + 
Vspirah turned on from the start, i.e., without some initial tran- 



sient time of growth of the non-axisymmetric part, and for just 
one value of the Jacobi constant. Figure llf shows the full 're- 
sponse model' obtained precisely via the method suggested by 
Patsis (2006). We thus take the initial conditions of 10^ particles 
on circular orbits of the axisymmetric part of the potential and 
on part of a uniformly populated disk, with distances in a range 
corresponding to the whole range of Jacobi constant values up 
to Ej^i4 of the full potential. The non-axisymmetric part of the 
potential was now introduced softly (in two pattern periods), and 
the particles were found to settle to a nearly invariant distribution 
in the configuration space after about 15 pattern periods. Figure 
llf shows this distribution at a time corresponding to the 25th 
pattern period. Clearly, the particles have settled to orbits sup- 
porting both the bar and the imposed spiral arms. But the most 
interesting feature of this distribution is that the maxima of the 
density of the response spirals depart from the maxima of the 
imposed spiral at nearly the same points (points A and A') where 
this happens for the invariant manifolds of Fig .4b, which are also 
essentially traced by the points in Figs. 1 lb,d, or e. The interest of 
this result Ues in that the initial conditions of the response model 
of Fig. llf are selected from a nearly homogeneous distribution 
in space, and thus they are by no means associated with the in- 
variant manifolds. In this respect, the invariant manifolds play 
for chaotic orbits a role similar to that of stable resonances for 
regular orbits, i.e. the manifolds are able to capture the chaotic 
orbits in their neighborhood and to create response spiral arms. 

4. Conclusions 

In the present paper we examined the applicability of the the- 
ory of the invariant manifolds emanating from the unstable 
short period orbits around the Lagrangian points Li and L2 of 
a barred galaxy in the self-consistent models of Kaufmann & 
Contopoulos 1996) which are rough models representing some 
features of real barred-spiral galaxies. We also tested the theory 
in weak or mean spiral versions of these models. Our conclu- 
sions are the following: 

1) When both the bar and spiral components of the self- 
consistent models are taken into account, the projection of the 
invariant manifolds on the configuration space produces a pat- 
tern that follows closely the imposed spiral pattern of the self- 
consistent model. 

2) The addition of the spiral potential produces two effects: 
a) it shifts the position of the unstable Lagrangian points Lj, or 
L2, both radially and angularly, and b) it enhances chaos locally, 
in a zone around corotation, due to the increase of the ampli- 
tude of the non-axisymmetric perturbation. The azimuthal dis- 
placement of the Lagrangian points is the most important factor 
for the manifolds to obtain a spiral form. In artificial models 
in which the whole non-axisymmetric perturbation is 'aligned' 
to the bar (e.g. by setting the phase of the spiral term equal to 
fl) = ff), the manifolds yield ring rather than spiral structures. If 
the amplitude of the perturbation is pushed to very high values, 
the lowest order lobes of the invariant manifolds can determine 
spiral patterns. However, we find that such patterns i) do not fit 
the spiral patterns of the self-consistent models, and ii) disap- 
pear after a few iterations of the calculation of the manifolds. 
This could probably be improved also by considering a lower 
n-value of the Ferrers bar, yielding a more gradual decrease of 
the bar force beyond the bar's end, and/or rendering the bar's 
outline more rectangular. 

3) We construct a simple theoretical model yielding the 
boundaries of the apocenters of the chaotic orbits in the con- 
figuration space. Since the manifolds are developed within these 
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boundaries, we demonstrate that they necessarily take the form 
of spiral arms, provided that the azimuths of the Lagrangian 
points Li, L2 are shifted with respect to the bar's major axis. The 
outermost radial limits of the manifolds are posed by rotational 
KAM curves outside corotation. However, in practice there are 
more stringent limits posed by the existence of cantori limiting 
the chaotic diffusion. The cantori define a narrow stickiness zone 
beyond corotation. The stickiness causes confinement of the in- 
variant manifolds inside this zone for times comparable to the 
age of the galaxy. 

4) A simple iterative method is proposed to calculate 'spiral 
response' models on the basis of the theory of the invariant man- 
ifolds, starting from the invariant manifolds of a pure bar model. 
We also discuss the relevance of the response models of Patsis 
(2006) to the theory of the invariant manifolds. Following a gen- 
eral property of the Hamiltonian flow in a connected chaotic do- 
main, a small sub-domain of initial conditions embedded within 
the chaotic domain yields successive Poincare maps following 
closely the form of the invariant manifolds of the main fami- 
lies of unstable periodic orbits in the same domain. The conse- 
quences of this effect are tested by calculating a full response 
model via the method proposed by Patsis (2006). Despite the 
fact that the initial conditions of the particles in the latter model 
have no relevance to the invariant manifolds, the maxima of the 
response spiral arms in the final state are found to follow closely 
the underlying pattern formed by the invariant manifolds. 
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APPENDIX 

A. Permissible regions for motion 

The Hamiltonian ([T]i, where Pr - r and pg - defines 
a curve of zero velocity r - -Q given by 



1q}/ + Vf,(r) + V,{r,G). 



(Al) 



where Voir) and V\{r,ff) are the total axisymmetric and non- 
axisymmetric potential terms respectively. Denoting by the 
corotation distance when only axisymmetric potential terms are 
taken into account, and Vs = V{)(?".s), we have Q.^^ = V'^r,. The 
curves of zero velocity around L4 and L5 can be found by ex- 
panding the radius and the Jacobi constant with respect to the 

corotation values, r = rs + (r- rj), £7 = ~~T^ + ^.s + f^- We find: 
2h = -M{r - r,f + -^{r - r,f + ... + 2Vi(r, 0) (A2) 



Q2 

P 



4D? 



{ks denotes the epicyclic fre- 



where M 
quency). 

Consider a very simple form for V\, e.g. V\ - Acos(20). 
For all values of the Jacobi constant between Ej - Eji^ and 
Ej - E jj^^ we have /i - A < 0, thus, since r - r s is small, we find 
from (A2) that M > 0, or 2Q.pl > 1. Omitting terms of third 
order in (r - r^), we find from Eq.(A2) 
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■ - = +1 [2A cos(26') - 2h\IM 



1/2 



(A3) 



If -A < h < A the motions are outside two banana-like curves 
around L4 and L5. In particular, for r - rs we have, around L4, 



+ j cos {hi A), while for = we have 



2(A - h)IM 



1/2 



(A4) 



Equation (A4) gives the intersections of the banana-like curves 
of zero velocity on the axis connecting the center to the stable 
Lagrangian points. The forbidden regions of the CZV can be 
compared with the corresponding forbidden regions of the sur- 
face of section r = 0, which is given by Eq.(fT5]l. In the latter case 
we have p] = rV = QV^, OV^ - 2Q.Q.pr^ -n 2(yo + Vi - Ej) = 0. 
Developing r, Q and Vo in powers of r - and omitting terms of 
order higher than the second we find the boundary of the forbid- 
den region: 



■M + 



M 
4f22 



(r - r^y + 2A cos(26i) = 2h . (A5) 



This gives 



2A cos(26') - 2/1 



1/2 



2Q„ 



(A6) 



Since 2D.pl Ks > 1, this boundary is larger than the boundary 
(A3). 



